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Abstract. Mergers of compact objects may lead to differ- 
ent astrophysical phenomena: they may provide sources of 
observable gravitational radiation, and also may be con- 
nected with gamma-ray bursts. Estimate of the rate with 
which such mergers take place are based on assumptions 
about various parameters describing the binary evolution. 
The distribution of one of these parameters - the kick ve- 
locity a neutron star receives at birth will strongly influ- 
ence the number and orbital parameters of compact ob- 
jects binaries. We calculate these effect using population 
synthesis of binary stars and show that the expected com- 
pact object merger rate changes by a factor of 30 when 
the kick velocity varies within its current observational 
bounds. 
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1. Introduction 

The gravitational wave detectors LIGO and VIRGO will 
soon be operational. Their completion brings up the ques- 
tion of the possible sources of gravitational waves and also 
their brightness and observable rate. The most often con- 
sidered sources of gravitational waves are mergers of com- 
pact objects; neutron stars or black holes. A number of au- 
thors have already considered these phenomena, and sev- 
eral estimates have been published ([Nfarayan et al., 1991 



Phinncy, 1991^ Tutukov and Yungclson, 1993; Portegics 



Zwart and Spreeuw, 1996| ) 

There are two approaches leading to the calculation of 
the merger rate. The first method, let us call it the ex- 
perimental approach, is based on the observational fact 
that we do see three binary systems of neutron stars: 
PSR1913+16 ([Taylor and Weisberg, 1~989|), PSR1534+12 



flWolszczan, 199 ip , and PSR2303+46 ( [Taylor and Dewey 
1988]). Based on this number, and considering various 



observational selection effects, like e.g. the pulsar beam 
width, one can estimate the number of such systems in our 
Galaxy. Given the observed orbital parameters for these 
systems one can then calculate the lifetime of each one and 



thus obtain the expected merger rate in the Galaxy. This 
approach suffers from several weaknesses: it is based on 
observations of only three objects, so small number statis- 
tics strongly affects the results. Moreover, the estimate of 
the rate is based on theoretical assumptions regarding the 
selection effects. These assumptions may lead to system- 
atic errors of uncertain value. 

The other approach that we can use, (let us call it a 
theoretical one), is based on population synthesis of bi- 
nary systems. All compact objects that will finally merge, 
must have their origin as ordinary stars in binary systems. 
First studies in that field were performed for Monte Carl o 
simulations of radio pulsars by Dewey and Cordes (1987 ). 
Thus it seems that by modeling the evolution of stellar 
binaries and the statistical properties of large ensembles 
of binaries one can calculate the expected population of 
compact object binaries. Analyzing this population and 
then considering the systems that will have merged within 
the Hubble time one can also estimate the present com- 
pact object merger rate in the Galaxy. This calculation 
requires the supernova rate in the Galaxy, and also the 
fraction of all stars that are in binaries. This "theoreti- 
cal approach" also suffers from the fact that it is based 
on several assumptions. In order to create a population of 
binaries one requires distributions of initial parameters, 
like the mass of the primary star or the initial mass ratio 
in the binary. Some stages in the stellar evolution require 
parameterization as well, e.g. mass loss and angular mo- 
mentum loss through stellar wind, mass exchange in the 
common envelope phase etc. Another parameter that may 
strongly influence the number of compact object binaries 
as well as their lifetime is the value of the velocity kick 
that a neutron star receives as a result of supernova ex- 
plo sion. The importance of this paramet er has been shown 
by Portegics Zwart and Spreeuw (1996 ). 

The measurement of the distribution of kicks in su- 
pernova explosions is not easy. One approach that can be 
taken is to measure velocities of pulsars and use this dis- 
tribution as the distribution of supernova kicks. Here one 
has to take into account selection effects, like for example 
the fact the fast neutron stars may leave the Galaxy and 
not be visible as radio pulsars. A comprehensive study of 
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this and other selection effects ( Arzoumanian et al., 1997) 
shows that a large fraction of neutron stars has veloci- 



ties above 500 km s . Ibcn and Tutukov (1996) argue 



that the velocities of radio pulsars can be explained by 
just taking into account the recoil velocity due to mass 
loss in supernova explosions in binaries, and neglecting 
a possible " natal kick" ! An independent study of veloci- 
ties of young pulsar through measuring offsets from the 



We restrict the range of the masses to the interval between 
WMq to 50M Q as we want at least the more massive star 
to undergo a supernova explosion. 

We adopt the following distribution of the mass ratio 

q 



§(q) oc const 



(2) 



center s of superno va remnants (Frail et al., 1994; Lyne 



and Lc rimer, 1994 ) confirms that the distribution initial 



tail. On the other hand Blaauw and Ramachandran (1998 J 
argue that the observed properties of pulsars can be ex- 
plained in a model with a single unique value of the kick 
Vkick — 200 km s . A similar result has been obtained by 



There is some unc ertainty regarding $. Portegies Zwart 
and Verbunt (1996 ) use $ cx (l + q)~ 2 , while [Tutukov and 



velocit }cs of neutron stars may have a large high velocity Yungclson (1994|) use $ o c q. Yet another distribu tion is 



Lipunov et al. (1997 ) who show that the observed pop- 
ulation of neutron stars is consistent with the kick ve- 
locities in the range 150 - 200 km s _1 . The distribution 
of the kick velocities has been parameterized by differ- 



ent authors: Portegies Zwart and Spreeuw (1996) used a 
Gaussian with the width of 450 km s -1 , while |Cordes and 
ChernC'ff (1997) proposed a weighted sum two Gaussians: 
80 percent with the width 175 km s _1 and 20 percent with 
the width 700 km s _1 . Here we attempt to investigate sys- 
tematically the effect of the kick velocity distribution on 
the compact objects merger rate. In section 2 we describe 
the model of binary evolution used for simulating the pop- 
ulation of binaries, in section 3 we show the results of the 
simulation and calculate the compact object merger rate. 
Finally we discuss our results in section 4. 

2. Population and evolution of binary systems 

While the evolution of a single star is only a function of 
its mass and metallicity the evolution of a binary presents 
a more complicated problem. It depends on the mass of a 
more massive component, the mass ratio of the smaller 
to the larger star q, and on the initial parameters of 
their orbit: e - the eccentricity and a the semi-major axis 
of the orbit. Following the general approach in the field 
we assume that these parameters are independent, i.e. 
the probability density can be expressed as a product: 
p(M,q,e,a) = ^ , (A/)<j>(g)S(e) r(a). Our population syn - 
thesis code is mainly based on Bcthc and Brown (1998|) 



The evolution of a single star is described by Eggleton 



et al. (198£ ). In the following all masses are in the units 
of solar mass M© . 

2.1. Distribution of the initial parameters 

Beyond 10M© the distributio n of the masses of the pr i- 
mary stars that we use follows Bcthc and Brown (1998 ), 



*(M) oc M 



-1.5 



(1) 



used in the simulation by Pols and Marinus (1994 ), who 
use a distribution slightly peaked for small values of q and 
falling down when q goes to unity. 

Initial binary eccentricity e is assumed to have value 
between and 1 and its initial distribution is taken from 



( Duquennoy and Mayor, 1991 ) : 



H(e) = 2e. 



(3) 



The distribution of the initial semi-major axis a used 
in population synthesis codes is flat in the logarithm, i.e. 



r(a) oc 



1 



(4) 



Both stars are initially massive main sequence stars, with 
the radii of at least 4-5 R , so we take minimum value of 
a to be 10 Rq . However if the sum of radii of two zero age 
main sequence stars exceeds 10 Rq, we chose the minimal 
separation to be twice the radius of the primary compo- 
nent. This should give the stars the space to evolve with- 
out merging before they go through their main sequence 
lifetime. Nevertheless some systems might be born with 
high eccentricity, and thus be in contact and merge form- 
ing a single very massive star. At this point we stop our 
calculation as we are not interested in this type of mergers 
in our study. 

Although very wide binaries with periods up to 10 10 
days are observed (Duquennoy and Mayor, 1991) and 



models with initial separation up to 10 R© are used 
(Portegies Zwart and Verbunt, 1996), we have set the 



maximal separation to 10 R© as wider binaries are very 
unlikely to survive two supernova explosions and form a 
compact object system. Our models showed that setting 
the upper limit to 10 6 R© does not change the results 
noticeably. 

We assume that the distribution of the kick velocity a 
newly born neutron star receives in a supernova explosion 
is a three dimensional Gaussian and parameterize it by its 
width a v . We vary a v to determine the dependence of the 
merger rate and properties of the compact object binaries 
population on this parameter. Alternatively, we draw the 



Other luthors have used similar distributions, e.g. Porte- 



gies Zwart and Verbunt (1996 ) use the exponent of —1.7 



natal kicks from a weighted sum two Gaussians (Dewey 



and Cordes, 1987 ): 80 percent with the width 175 km s 
and 20 percent with the width 700 km s _1 . 
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2.2. Evolution of binaries 

In our calculations we neglect the effects of stellar wind 
on evolution of binaries. 

Tidal circularization takes place in binaries for which 
the size of any component (or both components) is large 
in relation to their separation( Zahn, 197§| ). This happens 
when the stars evolve to the giant stage, and/or if the 
initial binary separation is small. Wc follow the prescrip- 
tion proposed by Portegies Zwart and Vcrbunt (1996); the 
circularization takes place if the stellar radius of one com- 
ponent is larger then 0.2 of the periastron binary separa- 
tion. The orbital elements (a, e) change with conservation 
of angular momentum until the new binary separation is 
5 stellar radii of the component in which tidal effects take 
place, or the binary is totally circularized (e = 0). 

Binary stars may undergo mass transfer and common 
envelope phases in different evolutionary stages. Depend- 
ing on the physical conditions this will result in the change 
of mass of each star, mass loss from the system, and also 
in change of the size and shape of the orbit. We describe 
below three schematic evolutionary paths and accompa- 
nying mass transfer types used in the code. 



2.2.1. Evolution for small mass ratios, (q < 0.8897) 

The more massive star (primary) evolves faster in the bi- 
nary. We use the following approximate formula for the 
main sequence lifetime of a. star with mass Ml T MS = 



where a is binary separation, Mi is the mass of the giant 
losing material, M2 is the mass of its main sequence com- 
panion and indices i, f correspond respectively to the val- 
ues before and after the Roche lobe overflow. The parame- 
ter (3 describes the specific angular momentum of material 



leaving binary (Pols and Marinus, 1994). The above equa- 
tion assumes that (3 is constant with time. The value of (3 
is uncertain; typically the values such as (3 > 6 or (3 = 3 



are used, for discussion see Pols and Marinus (1994). Here 



all calculations are performed for (3 — 6. The initial pri- 
mary, now either a helium star (after the mass transfer) 
or a giant goes supernova, leaving behind a neutron star. 
The initial secondary evolves, becomes a giant and may 
begin to transfer material to the newly formed neutron 
star increasing its mass, and possibly turning it to a black 
hole. We assume that the maximum mass of a neutron 
star is 2.2M Q . 

In describing accretion on a neutron star we follow 
Bcthe and Brown (1998] ) , and we will call this regime type 
II mass transfer. A neutron star accretes matter while 
moving in orbit through the extended envelope of the gi- 
ant. We use the Bondi-Hoyle accretion rate 



M = irpvR 



(7) 



where p is the density in the vicinity of the compact object 
of mass Mi, v is its velocity, and R ac = 2GM\jv 2 is the 
accretion radius. Consideration of the energy loss equa- 
tions lead to the final mass of the compact object ( Bcthc 



20 X l l) 6 i W 111 W . - yrs Aftnr this time- thr. primary and Brown, 1998|) 



evolves to the giant stage which lasts about 20 per cent 
of its main sequence lifetime, and increases its size. For 
q < 0.8897 the secondary is still on main sequence at the 
end of primary giant stage. If the giant primary overfills its 
Roche lobe mass transfer to the main sequence companion 
(star 2, with mass MX) and mass loss takes place. We de- 
note this regime by type /mass transfer. The mass transfer 
continues until the giant is stripped of its envelope, and 
thus the final mass of the giant will be just that of the 
i ts helium core. We use t he approximation M{ — 0.3 M] 



f (MX 
Mi = 2.4 x —I 



1/6 



M 



1 • 



(8) 



We assume that the giant loses its envelope so its final 
mass is just that of the helium core M/ = 0.3 Mj. The 
orbital separation follows from the energy integration i.e. 
equations 5.18 through 5.23 in (Bethe and Brown, 199S): 



a f 



0.145 x ai 



Ml 
Ml 



(9) 



( [Bcthc and Brown, 1998| ). A part of the envelope mass 
0.7 Ml is lost from the system while another part is trans- 
ferred to the companion. The fraction of mass transferred 
to the companion is proportional to the square of the mass 



ratio ( Vrancken et al., 1991; Bethe and Brown, 1998 ), so 
the mass of the companion after the mass transfer is 



Mi =Ml(q + 0.7q'). 



(5) 



where q is the initial mass ratio. 

The orbit was already circularized by tidal forces, when 
the giant was filling its Roche lobe, and now its size 
changes due to mass transfer an d mass loss. We describe 
the change in semi- major axis by ( Pols and Marinus, 1994 ) 



m( m| 
m[ mJ 




(6) 



Finally the initial secondary undergoes a supernova ex- 
plosion, and provided that the system survives, we obtain 
a compact object binary, consisting of either two neutron 
stars or of a black hole and a neutron star. 



2.2.2. Evolution of two stars for intermediate mass ratio, 
(0.8897 < q < 0.95). 

When q > 0.8897 the secondary is already a giant when 
the primary explodes as a supernova. The upper limit is 
explained in the next subsection. 

For the case of small orbital separations the evolution 
begins in a similar way as for small q described above. 
The evolution begins in a similar way as in the case of 
small mass ratio, and the primary may transfer mass to 
the secondary as described by type I above. When the 
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secondary becomes a giant the primary is already a helium 
star, stripped of its envelope. The secondary may transfer 
mass to the primary and we approximate this process also 
by type I. The primary explodes as a supernova and forms 
a neutron star, which may accrete from the secondary, 
provided that it is still a giant. The mass transfer in this 
phase is described by type II. The secondary undergoes a 
supernova explosion and if the system survives we obtain 
a compact object binary. 

In the case of large orbital separations the first stage 
of mass transfer from the giant to the main sequence star 
(type I) does not occur. The stars might begin to interact 
only when both of them are already in the giant stage. 
In this case they loose the common envelope. To describe 
this stage of evolution we assume that the giants loose the 
entire envelopes and become helium stars of 30 per cent of 
the initial giant masses. We then approximate the change 
in the size of the orbit by: 



Of 
a, 



1 + 



(Mi 



M\)M\{\-q*) 



a CE M{M{ 



(10) 



As above a is binary separation, Mi, Mi are the masses 
of primary and secondary. The indices i , / correspond to 
the values before and after common envelope phase, re- 
spectively, and ac e describes the efficiency of the orbital 
energy expenditure in the dispersal of common envelope. 
Equation ( p~0| ) was adopted from Yungelson and Tutukov 
(1993[j]and describes the change of the orbital separation 
a following decrease of the orbital energy which was used 
to eject a common envelope. They have also showed that 
0.6 < cice < 1, and we have performed our calculations 
for different values: acE = 0.4,0.6,0.8,1.0. Our calcula- 
tions show that the results (the production rate of the 
compact object binaries, and their properties) are almost 
indistinguishable for these three values, so for simplicity 
we present our results for one value of acE = 0.8. 

Similarly to the case of small orbital separations if the 
secondary is still a giant it may transfer mass to the neu- 
tron star formed in the supernova explosion of the primary. 
This mass transfer is described as type II. 

2.2.3. Evolution of two stars of almost equal mass, 
(q > 0.95) 

Two stars almost simultaneously (within 10% of their 
main sequence lifetime) leave main sequence and become 
giants. If their radii are large enough in relation to the 
periastron orbital separation, the orbit is tidally circular- 
ized. If the binary separation is small enough the giants 
may finally overfill their Roche lobe and enter a common 
envelope phase, which we describe by equation ([!(]). The 
two stars go supernova one after another, and if the system 
survives, we obtain a compact object binary consisting of 
two neutron stars. 



M, 



(a) 13.34 



(b) 13.34 



(c) 4.00 



(d) 1.4 



(e) 1.4 



(f) 2.31 



(g) 2.31 




M 2 A e 



8.25 37.6 0.89 



8.25 7.69 0.00 



11.8 0.74 0.00 



11.8 0.94 0.76 



11.8 0.39 0.00 



3.55 0.01 0.00 



1.4 0.02 0.58 



0.166 xlO 9 yrs 



Fig. 1. An example evolutionary path leading to forma- 
tion of a black hole neutron star binary. Units of mass are 
solar masses, and units of distance are astronomical units: 
(a) initial binary; (b) after tidal circularization when star 1 
expanded to giant size; (c) system after first mass transfer, 
star 1 explodes as supernova; (d) system after first super- 
nova; (e) after tidal circularization when star 2 expanded 
to giant size; (f) system after second mass transfer, star 2 
explodes as supernova; (g) system after second supernova: 
black hole neutron star binary. 

After each process which changes the orbital param- 
eters and sizes of the components we check if the com- 
ponents sizes are not larger than the periastron binary 
separation. 

2.3. Supernova explosion 

The supernova explosions are very likely to have an impact 
on the binary systems. We assume that in each supernova 
explosion a neutron star with mass of 1.4Mq is formed 
and we neglect the interaction of the companion star with 
the envelope ejected in the supernova explosion. 

We draw a random time for a supernova explosion in 
the orbital motion of the binary. We find the relative po- 
sition and velocity of the two stars for this moment. The 
exploding star is then replaced by a neutron star and we 
add the kick velocity to its orbital velocity, while the ex- 
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Fig. 3. Possible results of the evolution of a binary with the initial orbital separation a = 20R© and eccentricity e = 0.5 
for different values of the initial primary mass M and the initial mass ratio q. The top left panel shows the systems 
that were disrupted in the first supernova explosion. The solid line corresponds to the initial secondary mass 10A/q. 
The top right panel shows the systems in which the neutron star, born in the first supernova explosion, merged with 
the companion. The bottom left panel shows the systems disrupted in the second supernova explosion. The crosses in 
the bottom right panel show the neutron star binaries and the filled circles are the black hole neutron star binaries. We 
present one thousand of each type of binaries. In this calculation we used the Cordes & Chernoff (1997) kick velocity 
distribution. 



panding supernova shell carries away a part of the total 
momentum of the binary. If the energy of such a system 
is less than zero the system is bound and we calculate 
the parameters of a new orbit of the binary and its new 
spatial velocity. For a detailed discussion of the effect of 



sudden mass loss and a random kick velocity on binary 
parameters see Hills (1983 ). 

After the first supernova explosion we verify if the neu- 
tron star does not collide and merge with the companion, 
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Fig. 4. Possible results of the evolution of a binary with the initial primary mass 20 M & and the mass ratio q = 0.5. 
The top left panel shows the systems disrupted in the first supernova explosion, the top right panel shows the systems 
in which the neutron star, born in the first supernova explosion, merged with the companion. The bottom left panel 
shows the systems disrupted in the second supernova explosion. The bottom right panel shows the black hole neutron 
star binaries, and there are no neutron star binaries formed for this mass ratio (q — 0.5). We present one thousand of 
each type of binaries. In this calculation we used the Cordes & Chernoff (1997) kick velocity distribution. 



i.e. if the radius of the companion is larger than the peri- 
astron binary separation. 

2.4- Energy loss through gravitational radiation 

Once a compact object binary is formed its orbit will 
evolve because of the gravitational wave energy loss. Or- 



bital energy loss through radiation of gravitational waves 
becomes important once a compact object binary is in a 
tight and/or highly eccentric orbit. The evolution of the 
semi-major axis a and eccentricity of the orbit e, in a bi- 
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M, 



(a) 10.7 



M 2 A e 



10.2 2.00 0.50 



t ft 



(b) 10.7 ^^k^fcf^ 10.2 1.50 0.00 



7iT 



(c) 3.22 '''••'ig^T) 3.07 0.03 0.00 



(d) 1.4 



(e) 1.4 



3.07 0.02 0.67 



1.4 0.04 0.98 



1.207 xlO 6 yrs 



Fig. 2. An example evolutionary path leading to forma- 
tion of a neutron star neutron star binary formation: (a) 
initial binary; (b) after tidal circularization when both 
stars expanded to giant sizes, system ejects common en- 
velope; (c) system after common envelope phase, star 1 
explodes as supernova; (d) system after first supernova, 
now star 2 explodes as supernova; (e) system after second 
supernova: neutron star neutron star binary. 



nary emitting gravitational waves have been calculated by 
Peters (1964 ). The lifetime of a compact object binary is 



t 



■merg 



5c 5 a 4 (l - e 2 ) 7 ' 2 
256G 3 MiM 2 (Afi + M 2 ) 



1 



73 
24' 



37 
96* 



(11 



where a is the semi-major axis of the orbit, e is its eccen- 
tricity, and Mi, M 2 are the masses of the compact objects. 

3. Results 

The code described above allows to generate populations 
of compact object binaries and trace their statistical prop- 
erties. In this paper we mainly concentrate on the depen- 
dence of these properties as a function of the kick velocity 
distribution. In Figures [l] and |^ we present two example 
evolutionary paths leading to formation of compact ob- 
ject binaries: black hole neutron star binary and double 
neutron star systems. 

Let us first consider different evolutionary paths of the 
binaries depending on their initial parameters. We first 



consider a model with the Cordes and Chcrnoff (1997) 



kick velocity distribution. Tracing the evolution now de- 
pends on four parameters: M (primary initial mass), q 
(initial mass ratio), a (initial semi- major axis), and e (ini- 
tial eccentricity). In order to visualize the evolutionary 



0.8 



0.6 



0.4 



(ns,ns) <H 




H-£-£ (bh,ns)<" 
(ns,ns) >H 



Fig. 5. Population of the compact object binaries in the 
four categories: double neutron star binaries that merge 
within the Hubble time (empty squares), and these that 
do not (empty circles), neutron star black hole binaries 
that merge within the Hubble time (filled squares), and 
these that do not (filled circles). The error bars represent 
the counting statistics of the simulation. 



effects we fix two of them and present the types of sys- 
tems obtained in the course of the binary evolution. In 
Figure |^ we fix the initial orbital separation a — 20 Rq 
and eccentricity e = 0.5. The graphs are empty in the 
lower left part for which M x q < 10M©. This is the re- 
' gion for which the secondary is not massive enough to 
undergo a supernova explosion. However, in some systems 
for which the primary mass is M > 20M Q , the mass of 
the secondary is increased by accretion when the primary 
goes through the giant phase. Most of the systems end 
up either by disruption in the first supernova explosion or 
by a spiral in of the neutron star to the secondary (top 
panels) The remaining systems may be disrupted in the 
second supernova explosion. The compact object binary 
population (bottom right panel) is bimodal. The neutron 
star neutron star binaries are formed from the systems 
with q very close to unity, while the black hole neutron 
star systems are primarily formed when q is intermedi- 
ate. Thus the initial distribution of the mass ratio q has 
a strong influence on the production rates of these two 
types of compact object binaries. The number of systems 
shown in each panel does not correspond to the actual pro- 
duction rates. We present one thousand systems in each 
panel. The actual calculation produces 51% systems with 
the secondary not massive enough to undergo a supernova 
explosion, 40% systems torn after the first explosion, 6.2% 
of systems merged after the first explosion, 2.7% systems 
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Fig. 6. Distribution of the initial orbital compact binary systems parameters in the space spanned by the period P and 
eccentricity e, for four different Gaussian distributions of the kick velocity <j„ = km s _1 - top left panel, 200 km s _1 
- top right panel, 400 km s _1 - bottom left panel, and 800 km s _1 - bottom right panel. Each panel contains 1000 
objects. 



torn after the second explosion, and 0.35% of compact 
object binaries. 

In Figure |] we present the results of the evolution of 
systems for which we fix the initial primary mass and the 
initial mass ratio, while varying the the initial orbital pa- 
rameters a and e. The systems with small orbital separa- 
tions and high eccentricity merge in the early phase of the 
evolution and are not considered in this paper. Disruption 
after the first supernova explosion may occur to all sys- 



tems. The population of systems that end up merging after 
the first explosion, are disrupted in the second explosion, 
or form a compact object binary originate from the same 
region of the parameter space. One should note that be- 
cause of the circularization of orbits already in the initial 
stages of the binary evolution the final population is not 
very sensitive to the initial eccentricity. On the other hand 
the distribution of the initial orbital separation is impor- 
tant, as the compact object binaries originate in systems 
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with relatively small orbital separations (see bottom right 
panel in Figure ||) . One should note that all the compact 
object binaries shown in Figure |4| are black hole neutron 
star binaries. Double neutron star systems are formed only 
when q is nearly unity in our simulations. As before the 
number of systems shown in each panel does not corre- 
spond to the actual production rates, and show one thou- 
sand systems in each panel. The actual calculation pro- 
duces 12% systems born in contact, 43% systems with the 
secondary not massive enough to undergo a supernova ex- 
plosion, 34% systems torn after the 1st explosion, 0.7% 
systems merged after the first explosion, 0.50% systems 
torn in the second explosion and only 0.15% mergers. 



We present the dependence of the production rates of 
different types of compact object binaries on the width of 
the kick velocity distribution a v in Figure §. The num- 
ber of double neutron star systems that merge within the 
Hubble time increases with the kick velocity, while the 
production rate of black hole neutron star systems be- 
comes smaller. These two rates are nearly equal when the 
kick velocity is roughly that given by Cordes and Chcrnoff 
(19970 



In Figure ^ we present the distributions of the orbital 
parameters of compact object binaries for a few represen- 
tative values of the width of the kick velocity distribution 
a v . Each panel contains one thousand compact object bi- 
naries. To understand these plots let us first consider the 
properties of the population of objects just before the sec- 
ond supernova explosion in the case a v = Okms -1 . Some 
of the systems are wide, with the eccentricity varying from 
zero to unity, however, most of the systems populate a re- 
gion with eccentricities above ~ 0.45. These systems orig- 
inate in binaries with small initial mass ratios. A char- 
acteristic property of this group is a correlation between 
eccentricity and period. Objects in this group originate in 
systems with small initial value of q. The masses of the 
primary compact objects in these systems have a narrow 
distribution, the mass of the secondary after the accretion 
phase weakly depends on the initial mass before accretion 
- see equation ||, and is typically M ' 2 ~ 2.3 Mq. The mass 
of the secondary just before the second supernova explo- 
sion, is the mass of the helium core of secondary star and is 
in the range 3.0 M Q < M' 2 < 5.5 M . The orbital param- 
eters of such system after an explosive mass loss (second 
supernova explosion) depend only on the total mass loss, 
see e.g. Pols and Marinus (1994). In our calculation the 
newly formed neutron star has a mass of 1.4 Mq, thus the 
relative mass loss x = (1.4 + M' 2 )/(M' 2 + M' 2 ) is in the 
range 0.47 < x < 0.69. Systems that loose more than half 
of the mass (x < 0.5) are unbound. Other systems become 
eccentric and the eccentricity is given by e' = (1 — x)/x. 
Hence the lowest eccentricity a system can have after the 
second supernova explosion in our simulations is e' = 0.44. 




10 3 I0 4 10 5 10 s 10' 10= 10 s 10 10 10" 10 12 10 13 10 14 
L m,r S « [years] 

Fig. 7. Cumulative distributions of the t m rg - the life- 
times of compact object binaries for four different kick 
velocities: the solid line for a v = km s^ 1 , the short 
dashed line for a v = 200 km s , the long dashed line 
for a v = 400 km s _1 , and the dash dotted line for 
a v = 800 km s" 1 . Each line was plotted from a distrbution 
of 1000 objects. 

The relation between the new orbital period P' and the 
new eccentricity e' for different relative mass loss x is 



P = P 



(1 



(1 - e') 3 / 2 



(12) 



where Pq is the orbital period before the explosion. This 
relation defines the curved shape of the distribution in the 
P, e diagram. The objects populating the right hand side 
of the P, e plot, originate from systems with intermediate 
or large q. The intermediate q objects went through ac- 
cretion onto the neutron star (regime II of a mass transfer 
described above) and therefore similar reasoning as above 
applies to them. However systems with high initial value 
of q results in wide binaries (periods larger than 100 days) 
with different eccentricities. The compact object binaries 
in Figure ^ (top left panel) with eccentricities below 0.45 
and orbital periods smaller then 10 3 days originate in bi- 
naries of intermediate and large initial mass ratio. 

In the case of non zero kick velocity the systems with 
high q are very unlikely to survive. In fact there are only 
two such systems on the plot for a v — 200km s _1 , and 
none for higher velocities. The escape velocity in long pe- 
riod systems is low, and if such systems existed after the 
first supernova event, they have are very likely disrupted 
in the second supernova explosion. 

As the kick velocity is increased the shape of the dis- 
tribution in the P, e diagram also changes. In the case 
a v = 200 km s _1 there is a small fraction of high ec- 
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centricity, low period systems. This region of the param- 
eter space fills up as the kick velocity increases, see lower 
panel in Figure o where the case u v = 400 km s _1 and 
a v = 800 km s are shown. This is due to the fact that 
with the increasing kick velocity the long period systems 
are easier to disrupt and the fraction of surviving short 
period systems becomes larger. We note that the distri- 
butions sho wn in Figure p| are similar to those obtained 
previo usly ( Portcgics Zwart and Spreeuw, 1996|; Tutukov 



and Yi ngclson, 1993| ). 

In Figure [7| we present the cumulative distributions of 
the lifetimes of compact object binaries for the same set 
of kick velocities as in Figure |[ In the case of no kick 
velocity a v — km s _1 the distribution is bimodal: the 
systems with small q merge typically within the Hubble 
time (which we take to be 15Myr). However the systems 
with higher q remain in wide orbits and their merger times 
exceed the Hubble time. With the increasing kick velocity 
only the tightly bound systems survive. The lifetime of 
a system scales with the fourth power of the semi-major 
axis, and therefore the median lifetime decreases with in- 
creasing kick velocity. In the case of a v — 200 km s _1 
the median lifetime is « 4 x 10 8 years, and it decreases 
by a factor of four when the kick velocity is doubled. For 
the highest kick velocity velocity a v — 800 km s" 1 the 
median lifetime is only « 3 x 10 years. One should note 
however, that the distribution of t mclgc is very skewed, 
and we present it on a logarithmic axis. Thus there is a 
long tail extending to about the Hubble time. Most of the 
mergers take place in the Hubble time and only a small 
fraction of the total population lasts longer. 

Finally in Figure |S| we present f merge - the fraction 
of all binaries with the primary star more massive than 
IOM0 that produce a pair of compact objects in a binary 
system. We calculate f merge for a large range of the kick 
velocity distribution widths. To calculate each point in 
Figure || we calculated one thousand of compact binaries. 
This fraction falls down very approximately exponentially 
with the increasing kick velocity. We have fitted a modified 
exponential to this dependence 

/merge oc exp(-4.21 - 8.51 x lO" 3 ^ + 2.6 x lO -6 ^) .(13) 

The fit is shown by the dashed line in Figure || and is 
accurate to about 6%. Note that equation ( |l3| ) can be 
used to determine the merger fraction for any kick velocity 
distribution that can be expressed as a linear combination 
of three dimensional Gaussians. Equation ([l3]) can be used 
together with Figure || to obtain the production rates of 
any type of compact object binaries as a function of the 
width of the kick velocity distribution. 




400 600 
cr v [km s' 1 ] 

Fig. 8. Fraction of binary compact objects that merge 
within the Hubble time (15Gyrs), as a function of the 
width of the kick velocity distribution. We also show the 
fit discussed in the text. Each calculation produced 1000 
merging binaries. 

rate as 0.02/$^ year -1 , and the binary fraction as 0.5/t>j ra , 
where /sn and fbin are factors of the order of unity. 

In the simplest case when we assume that the star 
forming process has been going at the same rate through- 
out the history of the Galaxy we obtain the compact ob- 
ject merger rate 



r = 0.0001 x f SN fun f, 



merge 



year 



(14) 



where fmerge is expressed in percents. Taking values of 
J "merge from our Figure |^, we may calculate number of 
expected compact objects gravitational merging events, 
which for, let say, a v = 100 km s _1 will be 0.0002 per 
year which yields approximately 1 event per 5000 years 
per Galaxy. 

Since the compact object merger rate is directly pro- 
portional to fmerge it also depends exponentially on the 
kick velocity distribution width! The assumption about 
the constant star formation throughout the history of the 
Galaxy is in fact not really crucial. As we have seen most 
of the mergers take place within a few times 10 8 years even 
for rather small velocity kicks, therefore in reality we are 
only concerned about the star forming history in the last 
10 9 years. The star forming rate in the Galaxy has most 



The actual compact object merger rate in the Galaxy 



can be calculated given the observed supernova rate, the 
fraction of stars that exist in binaries, and assuming some 
form of the star formation history. Assuming that there 
is one supernova explosion every fifty years in the Galaxy 
and that binary fraction is 50%, we denote the supernova 



probably been constant over the last 5 x 10 9 years (Miller 
and Scalo, 1979| ) 



The calculation of the detection rate in gravitational 



wave detectors (Abramovici et al., 1992) requires a num 



ber of additional assumptions, like for example the galaxy 
density in our loc al and far neighborhood etc., for a 
discussion see e.g. Phinney (1991); Chernoff and Finn 
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(1993). Regardless of the assumptions the detection rate 
is proportional to the compact object merger rate in the 
Galaxy, provided that stellar populations are similar in 
other galaxies. When calculating the expected rates one 
has to take into account the masses of the systems that 
merge. The volume of the space in a flux limited sample of 
events scales as oc M 5 / 2 ( [Tutukov and Yungclson, 1993Q . 
Thus mergers of heavy objects like a neutron star and a 
black hole, or a pair of black holes are visible in a larger 
volume and may yield a similar observational rate despite 
the fact that they are not as frequent as the double neu- 
tron star mergers. 

4. Conclusions 

We have modeled the evolution of binary systems using a 
Monte Carlo code. The results of the simulations are con- 



The rates and also perhaps measurements of the "chirp" 
masses, M. = /i^Ms, where [i and M are the reduced and 
total mass of binary system - Chernoff and Finn (1993), 



sistent with the results of other codes ( 


Portegies Zwart and 


Yungclson, 1998 


; Lipunov et al., 1997; Portegies Zwart and 


Verbunt, 1996). 


Using this code we find that the merger 



rate of compact object binaries and consequently the de- 
tection rate in gravitational wave detectors falls approxi- 
mately exponentially with the width of kick velocity distri- 
bution. While the code that we use is far from describing 
all the details of binary stellar evolution we must empha- 
size that it produces similar results to the ones obtained 
elsewhere and in this work we only concentrate on the 
relative scaling of the resultant merger rate with the kick 
velocity in a supernova explosion. 

The exact shape of the kick velocity distribution is 
very difficult to measure. |Cordes and Chernoff (1997 ) and 



Bcthc and Brown (1998 ) use a distribution which is a 
weighted sum of two Gaussian distributions: 80 percent 
with the width 175 km s -1 , 20 perc ent with 700 km s _1 ; 
Portegies Zwart and Spreeuw (1996| ) use a Gaussian with 
the width of 450 km s _1 . On the other hand Iben and Tu- 



tukov (1996) argue that no velocity kicks are required at 
all, however the lack of pulsars in wide binaries suggests 
that at least a small kick of a few tens of km s _1 must 



be present (Portegies Zwart ct al., 1997). Thus, the veloc 



ity kick determination remains uncertain. Consequently 
the detection rate estimates in gravitational wave detec- 
tors may be uncertain by this amount. Approximating the 
kick velocity distribution by a single Gaussian profile and 
changing its dispersion, we calculated the merger rate for 
a wide range of velocity kicks. Changing the width of the 
assumed profile within the values proposed by other au- 
thors, namely from cr m i n = 0.0 to <r max = 500 km s 
results in decrease of the merger rate by a factor of 30. 
The expected number of compact object mergers varies 
by more than an order of magnitude when the kick veloc- 
ity goes from 200km s _1 (the value preferred in population 
studies) to 500km s" 1 (the measured in the observed pop- 
ulation of pulsars). 

The measurements of gravitational wave signals may 
thus allow some determination of the kick velocities. 



and their distribution will pose yet another constraint on 
the stellar evolution. 
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